namespace hms::swe
{

template<typename MeshType>
template<BlockLocation loc>
void SolverThread<MeshType>::computeNewState(
	const Block<loc>& block,
	Array3XsMap& qMap,
	Ref<ArrayXXs> dBlockNew,
	Ref<ArrayXXs> vBlockNew
){

	block.forEachRow( [&](Index y){
		auto q { block.rowMap_mem(qMap, y) };
		auto d { block.template rowMap<1>(dBlockNew, y) };
		auto v { block.template rowMap<2>(vBlockNew, y) };

		auto dCond = q.template topRows<1>() < SolverType::minWaterDepth;

		v = dCond.template replicate<2,1>().select(
			0, q.template bottomRows<2>().rowwise() / q.template topRows<1>()
		);
		d = q.template topRows<1>();
	});

}

} // namespace hms::swe
